Discrete breathers on symmetry- determined invariant manifolds 
for scalar models on the plane square lattice 

G.S.Bezuglova, G.M. ChechinO and P.P.Goncharov 

Department of Physics, Southern Federal University, Russia 
(Dated: January 28, 2013) 

Abstract 

A group-theoretical approach for studying localized periodic and quasiperiodic vibrations in 2D 
and 3-D lattice dynamical models is developed. This approach is demonstrated for the scalar mod- 
els on the plane square lattice. The symmetry-determined invariant manifolds admitting existence 
of localized vibrations are found and some types of discrete breathers are constructed on these 
manifolds. A general method using the apparatus of matrix representations of symmetry groups 
to simplify the standard linear stability analysis is discussed. This method allows one to decom- 
pose the corresponding system of linear differential equations with time-dependent coefficients into 
a number of independent subsystems whose dimensions are less than the full dimension of the 
considered system. 
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Introduction 



The problem of energy localization in discrete nonlinear Hamiltonian systems attracts 
a considerable attention over the past few decades. In the framework of this problem a 
particular interest represent stationary discrete breathers (DBs) — spatially localized and 
time-periodic dynamical objects. The history of discovery and discussion of different discrete 
breather properties can be found in a number of detailed review papers l|-|5|. Various 
analytical, numerical and experimental methods were used to study breathers in a wide 
variety of discrete Hamiltonian systems. However, the most papers on discrete breathers 
deal with one- dimensional chains, and much smaller number of articles are devoted to study 



these dynamical objects in 2D and 3D periodic structures |6Hl2 



lj-[l9| 



In the present paper, we discuss discrete breathers in 2D square lattice with one degree 
of freedom per lattice site. Such dynamical models are called scalar. 

Different physical interpretation can be given to a scalar model. For example, it has 
been used for describing transversal mechanical vibrations of plane lattice in (see Fig. 
1), charge vibrations in an electrical network of nonlinear capacitors coupled to each other 
with linear inductors in (see Fig. 2), etc. 
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FIG. 1: Mechanical model. 



In vector models, more than one degree of freedom is associated with every lattice cite. 
For example, one can consider mechanical vibrations of the plane lattice with x- and y- 
displacements of the particles along its surface. 

In many papers DBs in 2D and 3D lattices are investigated by some approximate methods 
such as "Rotating Wave Approximation" (RWA) etc., and only several works are devoted to 
study these dynamical objects with the aid of numerically exact procedures. For example, 
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FIG. 2: Electrical model. 



in the paper |6|], stationary DBs were constructed for scalar models on square and hexagonal 
lattices with homogeneous interparticle potentials of different even degrees. 

In the present paper, we develop a general group-theoretical method to obtain invariant 
manifolds admitting existence of localized excitations. Then the search of DBs can be done 
in two steps. 

1. Singling out the above mentioned symmetry-determined invariant manifolds. To this 
end, we must know only the structure and the symmetry group of the considered physical 
system and we don't take into account any information about interparticle interactions. 

2. Constructing time-periodic dynamical regimes on each of these manifolds. Obviously 
such regimes represent stationary discrete breathers. Unlike the first step, we now need to 
know differential equations describing a given dynamical model. 

Below we discuss the above approach using as examples the following two scalar dynamical 
models on the square lattice: 

ftj + 7 • C" 1 = te+i,; - '/,,)"' ; - - <,, ,./r 1 • 

%j + 7 • Qi,j = Qi+1,3 + Qi-i,j + Qij+i + Qij-i - ( 2 ) 

For both models i = 1..M, j = 1..N and periodic boundary conditions are assumed. In 
equilibrium state all q^j are equal to zero. Despite the fact that the dynamical variable qtj 
associated with the site can be of different physical nature, we often refer to it as the 
particle (atomic) displacement. 

Eqs. ([I]) describe dynamical model on square lattice whose "particles" interact with the 
nearest neighbours and with a substrate by homogeneous intersite and on-site potentials of 
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the same degree m. Eqs. ([2]) determine the Klein-Gordon model, i.e. they describe an array 
of Duffing oscillators with linear interaction between the nearest neighbours. In both models, 
the strength of the on-site potential relative to the intersite potential is characterized by the 
coefficient 7, and we study stability of discrete breathers with respect to this parameter. 

Note that the model was used in [6( without the on-site potential (7 = 0). 

Why we consider these two models? The model ([1]) does not admit the harmonic ap- 
proximation for m > 2, and therefore it has no phonon spectrum. On the other hand, there 
exist exact solutions to this model representing Rosenberg nonlinear normal modes (NNMs) 



24| (see, also, |23j). Existence of such solutions is guaranteed by the possibility to separate 



space and time variables as a consequence of homogeneity of the potential. In the dynamical 
regime corresponding to a given Rosenberg NNM, variables turn out to be proportional to 
one and the same function of time f(t), while spatial profile can be determined from a sys- 
tem of nonlinear algebraic equations. Each localized NNM represents a stationary discrete 
breather. Moreover, linear stability of these modes can be analyzed much simpler than the 
investigation based on the well-known Floquet method which one has to exploit in general 
case. 

The second model possesses phonon spectrum and turns out to be more complicated 
for studying. Indeed, we must solve nonlinear differential equations to construct discrete 
breathers and need in the Floquet method to analyze their stability. 

It is essential that for both models we search breather solutions on the same symmetry- 
determined invariant manifolds which are singled out by group-theoretical methods. The 
above invariant manifolds can be used to construct DBs for any other scalar models on the 
square lattice. 

For vector models and/or for different lattices, one can find the corresponding symmetry- 
determined invariant manifolds using the same group-theoretical approach. 

Let us emphasize that nonlinear differential equations of 2D and 3D dynamical models 
possess, as a rule, many different solutions and a symmetry-related classification is very 
desirable for any procedure of their construction. 

Among solutions which can exist on the symmetry-determined invariant manifolds dis- 
cussed in the present paper may be not only DBs, but various localized quasiperiodic dy- 
namical objects (some discussion of such objects in one- dimensional lattices can be found in 



251]). The stability analysis of quasiperiodic regimes represents serious difficulties since one 
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cannot use the Floquet method in such cases. 



In 



341 ]. we have developed a group-theoretical method for splitting the linearized dynam- 



ical equations near a given dynamical regime (periodic or quasiperiodic) with an arbitrary 
symmetry group. In some cases, this method allows one to simplify considerably the stability 
studying of a given dynamical regime, since high- dimensional linearized system is decom- 
posed into a number of independent subsystems of linear equations with time-dependent 
coefficients whose dimensions can be sufficiently small. 

The present paper possesses the following structure. The group-theoretical method for 
finding symmetry-determined invariant manifolds is discussed in Sec. HI Discrete breathers 
in scalar model ([Q) with homogeneous potentials of m = 4, 6, 8 degrees are considered 
in Sec. [TTJ while their stability is analyzed in Sec. IIHI Construction and stability of DBs 
in the Klein-Gordon model 2 are discussed in Sec. IIVI The group-theoretical method for 
simplifying stability analysis of periodic and quasiperiodic regimes in dynamical systems with 
discrete symmetry is presented in Sec. |Vl The results of the paper are briefly summarized 
in Conclusion. 



I. SYMMETRY-DETERMINED INVARIANT MANIFOLDS 
A. General discussion 

All possible solutions to a given system of differential equations L can be classified by 
subgroups Gj of the "parent" symmetry group Go consisting of all transformations of dynam- 
ical variables which transform this system to its equivalent form. In other worlds, Go is the 
group of invariance of the considered dynamical model. In particular case of a Hamiltonian 
system, Go can be treated as a group of all transformations under which its Hamiltonian 
turns out to be invariant. 

If we take a certain solution to the dynamical equations L and begin to act on it by all 
transformations g G Go, then some of these transformations conserve the solution — they 
form a subgroup G G G . On the other hand, acting on the above solution by other g G G , 
i.e. by g G, we obtain the orbit of m exact solutions to which our original solution belongs. 
The integer number m is the index of the subgroup Gj in the group Go- 

Thus, for the symmetry classification of all possible dynamical regimes in the given model 
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with parent group Go, we must consider all its subgroups Gj G Go- Each Gj determines an 
invariant manifold in the phase space of the dynamical system: if we select initial conditions 
belonging to this manifold and begin to integrate the differential equations of our model, we 
obtain a dynamical regime which never leaves the above manifold during the time evolution. 
Invariant manifolds corresponding to subgroups of the symmetry group of a given model we 
call " symmetry-determined invariant manifolds" . 

The idea of classification of dynamical regimes by symmetry-determined invariant mani- 
folds was taken as a principle in the theory of bushes of nonlinear normal modes in physical 



systems with discrete symmetry 



20 



2 if ] (see also 22j). A given m-dimensional bush is 



described by m dynamical variables and represents an exact solution belonging to a certain 
symmetry-determined invariant manifold of the considered dynamical model. 
In kl 
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26-28 



34j . we have obtained low-dimensional bushes of extended vibrational 



modes for dynamical models on various space-periodic structures. For example, in 26] this 
has been done for all possible 3D crystal structures described by any of 230 space groups. 
We have studied existence and stability of one-dimensional and two-dimensional bushes of 



NNMs in the FPU chains in 



27 



The group-theoretical approach for finding bushes of NNMs and for analyzing their sta- 
bility properties was described in detail in [22j. However, up to the present time, we applied 
this approach to study bush existence and stability for the case of extended NNMs. For this 
purpose, we need to make use of the subgroups Gj of the parent space group Go which include 
some translational symmetry elements. The primitive cell of the resulting spatial pattern of 
the vibrational state is some times greater than the primitive cell of the equilibrium state. 

In the present paper, we consider stationary discrete breathers which represent localized 
vibrational regimes. As a consequence, they are associated with point subgroups of the 
corresponding parent space group. 



B. Choice of the parent symmetry group 

To begin the procedure of finding symmetry-determined invariant manifolds, we must 
choose a parent symmetry group Go which embraces all transformations of dynamical vari- 
ables under which the considered mathematical model is invariant. Taking into account that 
the stationary state represents one of the solutions to the model dynamical equations and 
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that, therefore, the symmetry of this state must be a certain subgroup of Go, we , firstly, 
consider the spatial transformations under which the square lattice is invariant. 

It is well known, that the symmetry of square lattice is characterized by the plane space 
group C\ v = P4mm (in Schoenflies and in the International notations, respectively) and we 
can assume G = C\ v . 

In general, the parent group Go can be wider than the symmetry group of the stationary 
state. For example, several stationary states with the same energy and symmetry can 
correspond to a given Hamiltonian system and, as a consequence, the symmetry elements 
which transform these stationary states into each other must also enter the group Go- 

Moreover, the parent group Go can contain some additional symmetry elements due to the 
specific type of interparticle interactions in the considered model. For both models (jTJ), (j2J), 
all forces acting on particles turn out to be odd functions of the atomic displacements qij (the 
corresponding potentials represent even functions of their arguments). As a consequence, 
changing signs of all dynamical variables, which we denote by symbol P, also represents a 
transformation under which models ([I]), (j2J) are invariant and P must be incorporated into 
the parent group Go- 

Therefore, in our case, the symmetry parent group Go can be written as the direct product 
of the space group C\ v and the group P): 

G = Cl®(E, P). (3) 

The last group consists of two elements, E and P, where E is identity element. 

For scalar models on square lattice which are not invariant with respect to the transfor- 
mation P, we can assume Go = C\ v . 

If one find any additional transformations which don't change a given mathematical 
model, they can also be incorporated into the parent symmetry group Go- It such a case, we 
will obtain more detailed symmetry-related classification of the dynamical regimes because 
some more subgroups correspond to the wider parent group. 

C. Wyckoff positions 

As was already discussed, for finding invariant manifolds associated with stationary dis- 
crete breathers, we must deal with point subgroups of the parent symmetry group Gq ((3]). 
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Such subgroups (they don't contain any translational symmetry elements), can be found for 
any of 230 space groups in the standard textbooks on crystallography. One can find there 
a list of so called Wyckoff positions (WPs) or regular point sets in the primitive crystal cell 
which correspond to point subgroups of a given space group. WPs determine points in the 
crystalline lattice with different local symmetry and, in any crystal, atoms are distributed 
among these points. 

However, WPs characterize the local symmetry of different points in crystal lattice in 
equilibrium state. On the other hand, for the vibrational regime corresponding to a given 
DB, we deal with atoms displaced from their equilibrium positions. As a consequence, for a 
fixed instant, we have an atomic pattern whose symmetry is characterized by a subgroup of 
the local symmetry group of the crystal equilibrium state. Thus, we must find all subgroups 
of the crystal local (point) symmetry groups. 



For the square lattice with plane group C\ v , there are six WPs (see, e.g., [32|) and only 



three of them must be considered for our purpose. In the International notation, they read: 

la, 4mm (0,0); (4) 



lb, 4mm (-,-); (5) 

2c, mm ^0, -\ . (6) 

Hereafter, we will refer to these Wyckoff positions as WP-1, WP-2, and WP-3, respectively. 
The following information is given for each WP in (j4])-([6]): the multiplicity (1; 1; 2) and the 
label (a; b; c), the local point group (4mm; 4mm; mm) and the localization of the fixed 
point of this group in square primitive cell with unit edge (0, 0); (j^, ; ^0, \ j . 

The local groups of other Wyckoff positions turn out to be subgroups of the local groups 
listed in (JI])-(j6]), and we will automatically take them into account when select all subgroups 
of the groups specified in (jl])-©. 

The WP-1 and WP-2 possess identical point symmetry groups (4mm), but these groups 
are located in different ways with respect to the square cell of the considered lattice: the 
fixed point of the group 4mm of WP-1 lies at the corner of the primitive cell, while the 
group 4mm of WP-2 located in the center of the cell (see, Fig. [3]). Because of this reason, 



S 



subgroups of corresponding point groups of WP-1 and WP-2 generate different invariant 
manifolds (see Table 1). 

WP-1 WP-2 WP-3 



FIG. 3: Wyckoff positions in the primitive cell of square lattice. 



The point group (mm) of WP-3 is located at the middle of the cell edge. 
Note that WP-1 and WP-2 consist of only one point per primitive cell (their multiplicity 
is unity), while WP-3 consists of two points per the cell (see, Fig. [3]). 

D. Construction of symmetry-determined invariant manifolds 



In general case, dynamical equations of a given scalar model on the square lattice can be 
written in the form 

Qmxn = F(Qmxn)- (7) 

Here matrix 



Q 



( qn to ••• qin 

Q2l 122 ■■■ Q2N 



MxN 



\ 



(8) 



\ QMl QM2 ■■■ Qmn J 

determines the set of all variables qij = qij(t) corresponding to the MxN fragment of the 
lattice, while matrix F represents r.h.s. of the differential equations of the considered model 
(examples of such models are given by Eqs. ([I]) and fl2J)). 

As was discussed in Sec. II C[ we must consider all subgroups of local symmetry groups 
of WP-1, WP-2 and WP-3 demanding that the matrix (pattern) Qmxn is invariant with 
respect to these subgroups. As a result, we will obtain all possible symmetry-determined 
invariant manifolds. 

The restriction on Qmxn to be invariant under the action of a symmetry group G C Go 
leads to certain relations between the dynamical variables qij(t): some of them turn out to 
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be equal or differ from each other only in sign. Below, presenting a symmetry determined 
invariant manifold, we denote equal dynamical variables by the same symbols a(t), b(t), 
c(t), d(t) ... and don't point out explicitly their time-dependence. 

In the present paper, we study only strongly localized discrete breathers. This allows us to 
deal with relatively small fragments QmxN of each invariant manifold. We choose concrete 
values of M and N from the condition that amplitudes of qij (t) for peripheral sites must be 
much smaller than those for the breather core. 

Let us demonstrate a procedure for finding 3x3 manifold Q^ 3 invariant with respect to 
the point group G = C iv . 

This group consists of eight elements: four of them (gx, g 2 , g%, 54) are rotations by 
the angles 0°, 90°, 180°, 270° around z axis orthogonal to the lattice plane, while other 
elements (g$, g$, g?, gs) are reflections in the mirror planes orthogonal to 3D vectors 
(1,0,0), (1,1,0), (0,1,0), (-1,1,0). This four reflections we denote by the symbols a x , a y , a xy , 
o-- xy , respectively. [Detailed description of the group C± v one can find below in Eqs. (l53~I) .] 

Above listed symmetry elements, acting on the manifold 



' ?11 912 <?13 * 



o (1) 

<<*3x3 



Q21 Q22 Q23 
Q31 Q32 Q33 



(9) 



transpose its elements = in a certain way. In particular, for symmetry element g2, 
representing rotation by 90° around Z axis passing through the center of the square cell, we 
have 



#2^3x3 



^ Ql3 Q23 Q33 ^ 
Ql2 Q22 Q32 



(10) 



\qu q2i q3i J 

while for g 5 , representing reflection a x , we obtain 



' qi3 Qi2 qu ^ 



95^3x3 



Q23 Q22 121 

\ Q33 Q32 q3i J 

Note that g 2 and g$ can be considered as generators of the group G = C± v : 
elements of this group can be expressed via different products of #2 and g§. 



(11) 



all other 
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As a consequence, to obtain manifold, invariant with respect to the group G — C^ V) we 
can require that this manifold is invariant relative to gi and g§ only: 



&QSs = QSs. (12) 

&QSs = QSs- (13) 
Then from Eqs. (fl2|l. (|9|) and (flQ|l, we obtain 



9l2 = 921 = 932 = 923 = b(t), 
911 = 931 = 933 = 913 = c(t), 

922 = a(t). 

Here we have introduced three new variables a(t), b(t) and c(t), instead of nine old variables 
qijit), i = 1..3, j = 1..3. 

Similarly, we must take into account invariance condition f|T3|) comparing Eqs. (Q and 
(|TT|) . As a result, we obtain the following form of the manifold Q^xs invariant with respect 
to the point group G = C^ v : 

' c b c « 



o (1) 

<<*3x3 

\^c b c J 



(14) 



6 a 6 

c 6 c 

Let us note that Eq. (1131) does not produce any additional restrictions on the dynamical 
variables as compared to Eq. ( |T2i) . i.e. invariant manifold ( IT^tl is actually determined only 
by Eq. f|T2|) . This means that above obtained manifold Q^ 3 turns out to be invariant not 
only relative to the group C4 = {C4} with one generator #2, representing rotation by 90°, 
but simultaneously relative to its supergroup C± v = {C4, a x } with two generators — g<i and 
g$. (Hereafter, we define any symmetry group by the list of its generators given in curie 
braces) . 

However, if we consider larger fragment of dynamical variables on square lattice, for ex- 
ample, the manifold Q^si Eq. (jTBl) indeed produce additional restrictions on some variables 
qij and, as a consequence, the manifolds invariant with respect to the groups C4 = {C4} and 
C^ v = {C^,a x } prove to be different (see, and in Table 1). 

Among point subgroups of the parent symmetry group, there are such that generate 
invariant manifolds admitting existence of localized, as well as delocalized dynamical objects. 
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Indeed, the manifold (fH|) allows localized periodic or quasiperiodic modes if \a\ » \b\ » 
\c\, while the manifold 

' b a 

b a c 
yb a c J 

does not admit localization in Y-direction because variables g^, expressed via a, b, c, have 
no tendency to decrease by amplitude from its center to the periphery. 

Note that some invariant manifolds can correspond to multibreathers of complex struc- 
ture, for example, with zero amplitude in the breather center. In Table 1, we present all 
invariant manifolds associated with subgroups of the parent symmetry group ([3]) which admit 
construction of simple discrete breathers — they are localized in both X- and Y- directions 
and have no zero in the center. 

However, at the end of each WP section of Table 1, we give the list of all other subgroups 
associated with this regular point set. Such information can be used for construction some 
different types of invariant manifolds, corresponding to this WP. 

Moreover, below each manifold from Table 1, we give its symmetry group which can be 
useful to expand a given manifold in the case of bad discrete breather localization (vibrational 
amplitudes slightly decrease from the breather core to its periphery). 

In Table 1, we present only small fragments of invariant manifolds QV' (j = 1..27), for 
example, for manifolds corresponding to WP-1 we show 5x5 fragments and give their 3x3 
parts in frames. 

Note that manifolds with a minus sign before some variables correspond to subgroups 
whose elements contain operator P (independently or in combination with other symmetry 
elements). We list these manifolds after those corresponding to the parent group G = C± v . 
The former can be used only for models with even potential, while the latter correspond to 
the models with an arbitrary potential. 

Finally, let us emphasize once more that all invariant manifolds presented in Table 1 
admit existence of periodic, as well as quasiperiodic dynamical objects. 

In the following sections of the paper, we discuss construction of DBs as spatially localized 
and time-periodic vibrations on some manifolds listed in Table 1 for the dynamical models 
© and P). 
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Table 1 

Symmetry-determined invariant manifolds for scalar models on the square lattice. 

WP-1 



e 


d e 


c 


b c 


b 


a b 


c 


b c 


e 


d e 


) : 


{C 4 , 


i 


9 3 


d 


b e 


b 


a c 


e 


c f 



\ 



y p k h I m j 
Q (6) : Wxy} 



/ g e d f g\ 




\g f d e g J 
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II. CONSTRUCTION OF DISCRETE BREATHERS FOR THE MODEL WITH 
HOMOGENEOUS POTENTIAL 

Each invariant manifold from Table 1 depends on a number of arbitrary parameters (a, 
b, c, ...). To construct a discrete breather, we must find such values of these parameters 
which lead to a spatially localized and time-periodic vibration when they are used as initial 
values for integrating dynamical equations of a given mathematical model. 

In this section, we discuss construction of DBs for the model (CQ) which corresponds to 
the case of homogeneous potential of m-degree. Similar model without on-site potential and 
without studying breather stability was considered in jf5|. Unlike this paper, we study the 
dependence of existence and stability of DBs with respect to the relative strength (7) of 
on-site and intersite parts of the potential energy of the model ([1]). 

The specific structure of Eqs. (CQ) admits the space-time separation and, as a consequence, 
we can treat discrete breathers for the case of the homogeneous potential of m degree in 
terms of localized nonlinear normal modes (NNMs) introduced by Rosenberg 24J. To this 
end, we assume that 

qi)j (t) = k ir f{t), (15) 

where i = 1..M, j = 1..N, while k^ are constant coefficients. 

Substituting the ansatz ( IT5l) into differential equations (CQ) and requiring that all they are 
equivalent to each other, we obtain a number of nonlinear algebraic equations, which deter- 



14 



mine the spatial profile of DB, and one ("governing") differential equation, which determines 
time-dependence of all the dynamical variables qij(t). This time-dependence is described by 
the single function f(t). 

If we now take into account particular structures of the invariant manifolds from Table 1, 
the number of unknown coefficients kij will be equal to the number of manifold arbitrary 
parameters a, b, c, ... minus one, since one of these parameters can be fixed. 

In {25), Q' we have used Rosenberg modes technique for constructing DBs and studying 
their stability in one-dimensional lattices. In the present paper, we use the same technique 
to investigate DBs in the scalar model ([I]) on the two-dimensional square lattice. 

Let us illustrate the procedure of the Rosenberg mode construction using as an example 
Q3X3 invariant manifold (|14j) . Only three independent dynamical variables, a(t), b(t) and 
c(t), correspond to this manifold. According to the definition of Rosenberg modes, we 
assume that these variables are proportional to the same time-dependent function f(t): 

a(t) = a-f(t), b(t) = b-f(t), c(t) = c-f(t). 
Here a, b, c in r.h.s. of the equalities are constants which we use, for simplicity, instead of 
the coefficients from Eq. ( fl5l) . Therefore, the Rosenberg mode R(t) can be written in 
the form: 



R(t) 



' c b c ^ 

b a b 
y c b c j 



/(*). ( 16 ) 



where the vibrational amplitude corresponding to the center site of the manifold can be 

assumed equal to unity (a = 1). Substituting dynamical variables %■(£) from (fT6|) into 
differential equations ([1]), we can obtain the following algebraic equations 

b ■ [-7 + 4(6 - l)™- 1 ] = -7 • b m - x + 2 • (c - b) m - x + (1 - b) m -\ (17) 

c • [-7 + 4(6 - l) m - 1 } = -7 • c m ~ l + 2 ■ {b - c) m -\ (18) 
while the governing equation reads 

/(*) + P 2 ■ =0, P 2 = 7 - 4 • (b - a) m -\ (19) 
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As a rule, there are many different solutions to Eqs. (ITT]) . (ITSjl and we must demand that 
\a\ > \b\ > |c| to provide a localized breather profile. 

There are many methods for solving nonlinear algebraic equations, such as various ver- 
sions of the Newton-Rafson method, the steepest descent method, etc [5]. We use for this 
purpose the standard procedure of the mathematical package MAPLE. In this way, we have 
obtained breathers on a number of invariant manifolds presented in Table 1 for different 
degrees m of the homogeneous potential and for a dense set of values 7 G [0; 10]. 

Let us consider some examples. 

1. For the homogeneous potential of m = 4 degree, we have obtained the following 
breather spatial profiles on the invariant manifold Q^ 3 for different values of 7: 

7 = 0: a = l, b = -0.25439, c = 0.00439; 

7 = 3: a = l, 6 = -0.17354, c = 0.00113; (20) 
7 = 6: o=l, b = -0.12261, c = 0.00032. 

We see that the coefficients b and c decrease with increasing 7, i.e. breather localization 
becomes better with increasing the strength of the on-site potential with respect to the 
intersite potential. (Note that the breather for 7 = coincides with that found in 6J). 

The time- dependence of the breather solutions with spatial profiles ([20]) determined by 
the governing equation 

f{t) + p 2 ■ f(t) 3 = 0, p 2 = 4(a - b) 3 + 7, (21) 

with analytical solution of the form 

/(t)=cn(w*,-j=), u = p-f(0). (22) 

In Fig. Hk, we show the coefficients b = 6(7) and c = 0(7) as functions of the parameter 
7 G [0;10]. Note that these functions possess opposite signs and are plotted on different 
scales. 

Thus, we can speak about the breather family because there exists a certain DB for 
any fixed value of the parameter 7. Moreover, one can see from Fig. Hk. that 6(7) and 0(7) 
represent, in the considered case, sufficiently simple functions which can be approximated by 
some polynomials with good accuracy. Using the mean-square method we find the following 
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FIG. 4: Discrete breathers for the model 1 on the invariant manifolds a) ^3x3! b) Q3X3, Q^x3' 
{7) 

^3x3- 



approximations: 

6( 7 ) = -0.133940 • 10- 2 7 2 + 0.297771 • 10 _1 7 - 0.252179, 

(23) 

c( 7 ) = -0.980247 ■ 10~ 5 7 3 + 0.213750 ■ 1(T 3 7 2 - 0.158541 ■ 10~ 2 7 + 0.423062 ■ 10~ 2 . 

The above discussed breather family on the manifold Q^3 f° r m = 4, 7 G [0; 10] is 
described by two variables 6(7), 0(7) because the third variable, a (see Table 1), can be 
considered as a constant: a = 1. Another fixed value of this variable leads to the breather 
spacial profile which turns out to be proportional to that for a = 1 (despite the algebraic 
equations (|T7|) . ( TT8|) for finding Rosenberg nonlinear normal mode profiles are nonlinear!). 

However, the time dependence of the breather solution, in particular, its frequency, does 
depend on the parameter a, since p = p(a) according to Eq. ( l2Tj) . 

Let us represent profiles for some other breathers whose stability we will discuss in the 
next section. 

a) Invariant manifold Q^x3 {CA,a y }. 



m = 6 : 0(7) = 1, 

6(7) = -0.733310 • 10- 3 7 2 + 0.228270 • 10^7 - 0.251758, 

c( 7 ) = -0.346820 • 10- 6 7 3 + 0.780072 ■ 10~ 5 7 2 - 0.593555 • 10" 4 7 + 0.157500 • 10~ 3 . 

(24) 



17 



m = 8 : 0(7) = 1, 
6( 7 ) = -0.132486 • 10~ 3 7 2 + 0.151641 • 10 _1 7 - 0.251854, 
c ( 7 ) = -0.103097 • 10- 7 7 3 + 0.255021 • 10~V - 0.215981 • 10~ 5 7 + 0.636736 • 10~ 5 . 

(25) 

b) Invariant manifold Q^ 3 {C 4 ,ay}. 

m = 4 : 0(7) = 1, 

6(7) = 0.206049 • lO^V - 0.462088 • 10~ V + 0.864977 • 10 _1 7 - 0.539657, 

c( 7 ) = -0.155129 • 10- 3 7 3 + 0.332604 • 10~V - 0.232426 • 10^7 + 0.538939 • 10" 1 

(26) 

m = 6 : 0(7) = 1, 

6( 7 ) = -0.334370 • 10~ 3 7 3 + 0.420770 • lO'V + 0.272984 • 10^7 - 0.500903 

c( 7 ) = -0.541969 • 1(T 5 7 3 + 0.153732 • 10~ 3 7 2 - 0.143670 • 10~ 2 7 + 0.445843 • 10~ 2 

(27) 

m = 8 : 0(7) = 1, 

6( 7 ) = 0.104004 • 10" 3 7 3 + 0.174863 • 10~ 3 7 2 + 0.151320 • 10 _1 7 - 0.500456, 

c ( 7 ) = 0.161985 • 10~ 6 7 3 + 0.133814 • 10~ 5 7 2 - 0.751460 • 10~ 4 7 + 0.460920 • 10~ 3 

(28) 

b) Invariant manifold Q^l {C^A "?/}- 

m = 4 : 0(7) = 1, 

6( 7 ) = 0.479180 • 10- 4 7 3 - 0.131456 • lO^V + 0.161263 • 10 _1 7 - 0.1356280688, 

c ( 7 ) = -0.971222 • 10~ 6 7 3 + 0.208977 • 10~ 4 7 2 - 0.155620 • 10~ 3 7 + 0.442289 • 10~ 3 . 

(29) 

m = 6 : 0(7) = 1, 

6(7) = -0.259370 • 10- 4 7 2 + 0.113473 • 10" 2 7 - 0.344276 • 10" 1 , (30) 

0(7) ~ 10- 7 . 
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m = 8 : a(-y) = 1, 

b(j) = -0.479130 ■ 10~V + 0.659596 • 10~ 4 7 - 0.813363 ■ HT 2 , (31) 

0(7) ~ 10~ 16 . 

In Fig. HJd, one can see an example of nontrivial profiles of DBs for the homogeneous 

potential of m = 8 degree. These breathers are constructed on the invariant manifold 

(7) 

Q$x3 which is characterized by five arbitrary parameters b, c, d, e, f (the parameter a = 1 
corresponds to the breather centre). Parameters e and / corresponding to the lattice sites 
more distant from the central site with a = 1 are very small (e, / ~ 10 -6 ) as compared to 
parameters b, c, d associated with the nearest sites. 

There is a critical point 7 = 7 C ~ 5.885 in which some DBs of different symmetry 
coincide with each other. Indeed, the dotted line on which b = c = d corresponds to the 
breathers on the invariant manifolds Q!ix3 w hh the symmetry group {C4, o~ y } = C^ v . The 
solid lines (c = d 7^ b) corresponds to breathers on the manifold with the symmetry 

group {C 2 , <J y } = C 2v , while the dash-dot lines (b = d 7^ c) corresponds to DBs on the 

(7) 

manifold Q 3 ^ 3 whose symmetry is described by the group {cr y } = C s . 

Note that the above mentioned breathers, except for those who possess the symmetry 
group {C4, a y } = C iv , turn out to be unstable (see below). 

III. STABILITY ANALYSIS OF DISCRETE BREATHERS FOR THE MODEL 
WITH HOMOGENEOUS POTENTIAL 



Discrete breathers represent time-periodic dynamical regimes and, therefore, standard 
Floquet method can be used for analyzing their stability. However, for the scalar model 
(CQ) on the square lattice with m-degree homogeneous potential the stability studying of 
discrete breathers can be essentially simplified, if they are Rosenberg nonlinear normal 
modes (in the next section we will consider breathers in the model 1 which are not such 



modes!). In 25j, |30j, we have demonstrated this idea for the case of one-dimensional chains 



with m = 4 homogeneous potential, while now we use the same method for analyzing the 
breather stability in the two-dimensional scalar model ([1]) described by the homogeneous 
potential of arbitrary (even) degree m. 

Let us consider a localized nonlinear normal mode. Being periodic in time, this dynamical 
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object represents a certain discrete breather. Let vector 8(t) determines a set of infinites- 
imal deviations of all the particles from the exact breather solution R(t) [for example, see 
Eq. ( 1201) ]. Linearizing nonlinear differential equations of the model (CQ) with respect to all 
components of the vector 6, we obtain the following system of linear differential equations 
with time-periodic coefficients (some computational details can be found in [251]]): 

5 = (m - 1) ■ f rn - 2 B-8. (32) 

Here B is a constant symmetric matrix whose dimension, as well as the dimension of the 
vector 8, is equal to the number of all sites of a chosen lattice fragment, while f(t) is the 
solution to the governing equation for the Rosenberg NNM. 

For the model ((Tj), described by the homogeneous potential of m degree, the time-periodic 
function f(t) is the solution to the equation 

m+p'-fitr-'^o, (33) 

with initial conditions 

/(0) = A, /(0) = 0. (34) 

Here, A is the breather amplitude (displacement of the central particle), while p 2 depends 
on the breather profile. For example, for a breather on the invariant manifold Q;sx3 ^ n the 
case of the model 1 with m = 4, we have Eq. (J21]). 

In the above approach, investigation of stability of the considered dynamical regime 
(discrete breather, in our case) is reduced to analyzing the stability of the zero solution of 
the linearized variational system f[3"2"]) . 

Note that the system (132]) possesses a very specific structure: it contains the constant 
matrix B with time-periodic coefficient (m — 1) • f m ~ 2 (t) standing in front of it. 

Now we can reduce the symmetric matrix B to a diagonal form D with the aid of a 
certain orthogonal transformation 

S T BS = D. (35) 

(S T is a transpose of S). 

Eigenvalues Xj of the matrix B are diagonal elements of the matrix D, while eigenvectors 
of B form the columns of the matrix S. 

Introducing a new infinitesimal vector z instead of the old one 6 with the aid of the 
transformation 6 = Sz, we decompose the system (132]) of coupled equations into independent 
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scalar equations of one and the same type: 

z=(m-l)'X r fit)" 1 ' 2 ■ z. (36) 

With the help of scaling 

f{t) = Ai,(t), t = -^ (37) 
one can rewrite Eqs. f l33]) . (JM]) and fl36|) in the form 

< + ^ m - 1 (r) = 0, (38) 

#)) = 1, Vr(0) = 0, (39) 

4' = ^ m - 2 (r)-A r z, (40) 

z(0) = 0, 4(0) = 0, (41) 

where 

A, = (42) 

Here we denote differentiation by the new time variable r by prime. Eqs. (j4ip appear 
because we study stability of the zero solution of the differential equation (j4"0~|) for the 
function z(t). 

We will refer to the parameters Aj as stability indicators. Breather solution proves to be 
stable if all Aj fall into intervals of stability of the zero solution to Eq. ( 1401 . 

Note that one of the stability indicators, say A , always belongs to the boundary between 
the first region of instability and the second region of stability of zero solution to Eq. (1401) . 
The infinitesimal vector corresponding to this indicator turns out to be exactly proportional 
to the breather profile (i.e., it is located along this profile) and, therefore, such indicator A 
does not affect the breather stability. 

For m = 4, stability analysis of zero solution of Eq. (1401 can be fulfilled with the aid 



of the Lame equation (see, 2510 . For m > 4, such analysis represents essentially more 



complicated analytical problem 35(. However, some numerical methods can be efficiently 
used for studying stability of the zero solution of Eq. (140|) . Using these methods we were 
able to reveal that discrete breathers prove to be linear stable if all indicators Aj fall into 
the following intervals: 

[0; 1], [m-l;m + 2], [3m - 2; 3m + 3], ... (43) 
21 



Let us consider stability of discrete breathers listed in Eqs. (j23]) - (l3Tj) . The indicators 
Aj = are certain functions of the parameter 7 describing relative strength of the on- 

site and intersite potentials of the model 1. In Fig.0 we present functions Aj^), 7 G [0; 10] 
for different manifolds and different degrees m of the homogeneous potential. 



Q (1) {C 4 , aj 



A. 




Q (8) {C 4 , a y } 

A. 



m=4 




Q (15) {C 2 P, a y } 



A. 



A 



A. 




m=4 



m = 6 



FIG. 5: Stability indicators Aj(7) for discrete breathers on the invariant manifolds {C^,a y }, 



{C 4 ,a y }, Q (lb ) {C 2 P,a y } for the model 1. 



1(15) 



Note that stable and unstable regions of the zero solution of Eq. fHOj) alternate according 
to Eq. (I43p . The first stable region [0, 1] are one and the same for all m, while the width 
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of the first unstable region [l,m — 1] increases with increasing m. Because of this reason, 
we depict in Fig. [5] only lower part of the first region of instability (it is presented in grey 
color) . 

Let us comment in these figures. 

1) DBs corresponding to the invariant manifold Q^ 3 with the symmetry group {C4, a y } = 
C4V centered at the corner of the primitive cell (WP-1) are described for m = 4, 6, 8 by Eqs. 
( I23TI25"]) . Nine degree of freedom are associated with this 3x3 fragment of the square lattice 
and, therefore, there are nine indicators Aj which determine stability of breathers on Q 3 ^ 3 - 

For example, for m = 4, 7 = 0, we have found the following indicators Aj (j — 0, 1, 8): 

A = 3, Ai = A 2 = 0.65095, A 3 = 0.64879 

(44) 

A 4 = 0.09128, A 5 = 0.05089, A 6 = A 7 = 0.04873, A 8 ~ lO" 10 

Note that some indicators Aj are equal to each other (we consider the source of such 
degeneration elsewhere). The functions Aj = Aj(7) are depicted for the breathers on Q^ 3 
for m = 4,6,8 7 G [0; 10] in Fig HO For m = 4, all stability indicators Aj(7) fall into the 
first stability region since < Aj(7) < 1, j — 0..8, 7 G [0; 10]. Note that we cannot see all 
these indicators in Fig. [5] because some of them are too small. Therefore, discrete breathers 
described by Eq. fl23|) (m = 4) are stable dynamical objects for all values of the parameter 
7 on the considered interval. 

However, one can see essentially different stability properties of discrete breathers for 
m = 8 on Q3X3. Indeed, one curve, Ai(7), lies out of the first stability region for 7 < 7 C « 
5.8848. This means that DBs for small strength of the on-site potential turn out to be 
unstable. On the other hand, for 7 > j c , this curve enters the first stability region and DBs 
become stable for 7 greater than j c . 

Finally, for m = 6, we have an intermediate case between two above discussed cases: the 
curve Ax (7) begins exactly from the upper end of the first stability region and, therefore, 
the breather possesses marginal stability for 7 = 0. All DBs for 7 > are stable since 
A,( 7 ) < 1, j = 0..8, 7 e(0;10]. 

2) Another situation takes place for DBs (1261- 1281) on the invariant manifold Q^s with 
the group {C 4 ; a y } = C 4 t> corresponding to the center of the primitive cell (WP-2). Indeed, 
we see in Fig. [5] that for m = 4, m = 6 and m = 8, some stability indicators lie outside 
the first stability region up to certain critical values 7 C of the parameter 7. These critical 
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values 7 C = 7 c (m) are different for different m. Moreover, we see that 7 c (m) increases with 
increasing m. In other words, we need more strong on-site potential for stabilizing DB for 
the greater m. 

3) Finally, for DBs determined by Eqs. f r2"9~j) -( l3Tj) on the invariant manifold Q^l with the 
symmetry group {C2P; o~ y } = C2V corresponding to WP-2, we find Aj(-f), j = 0..15 which 
also are depicted in Fig. [5j From these plots, one can see that DBs for m = 4 are unstable 
for all values of 7 G [0; 10]. On the other hand, for m = 6 and m — 8, DBs turn out to 
be stable for small values of 7 up to the critical values 7c ~ 5.3104 and ~ 19.1248, 
respectively. In all above discussed cases, there is the same tendency: a sufficiently strong 
on-site potential can stabilize breather vibrations. 

In conclusion, let us pay attention to one important property of the discrete breather 
stability for the models with homogeneous potential of an arbitrary degree m. In such 
models, the stability of the breathes which represent Rosenberg nonlinear normal modes 
don't depend on their amplitude. This fact becomes obvious, if we take into account that 
the amplitude A does not enter Eqs. (I38l)-( [4*2l) which fully determine the stability of the 
zero solution to decoupled linear equations derived from the system ( 1321) . 



IV. DISCRETE BREATHERS FOR ARBITRARY SCALAR MODELS ON THE 
SQUARE LATTICE 

A. Construction of DBs for the model 2 

The scalar model 1 with homogeneous potential of m > 2 degree discussed in the previous 
section represents a sufficiently exotic case from the physical point of view. Indeed, there 
is no phonon spectrum in this model since it does not admit the harmonic approximation: 
Taylor-series expansions of r.h.s. of Eqs. ([Q) with respect to all dynamical variables 
don't contain any linear terms. On the other hand, in the more realistic models, admitting 
the harmonic approximation, there are no localized Rosenberg nonlinear normal modes. 
Therefore, we have to use different and more complicated methods, as compared to those 
described in Sec. IHIl for constructing DBs and for stud ying their stability. The most widely 
known method for constructing DBs was developed in [29| . It based on the Newton- Rafson 
iterative scheme for finding initial values of all dynamical variables which lead to localized 
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and time-periodic solution of the dynamical equations describing the considered model. 
However, such approach need, as a rule, very good initial approximation for the breathers 
profile at t = and, even in this case, some problems with convergence of the numeric 
procedure in the many-dimensional phase space PS can arise. We prefer to construct DBs 
with the aid of a certain variant of the steepest descent method in the above mentioned space 
PS (see, for example, j^]). Using this method, we minimize the sum of square deviations 
between all dynamical variables at the instants t = and t = T where T is an arbitrary fixed 
oeriod of the breathers solution. In some cases, we used the "pair synchronization method" 



30t ] in a manual regime for constructing plausible initial conditions. Runge-Kutta methods 
(rkf45 and dverk78) were applied for integrating differential equations of the considered 
model from t — to t — T. Such approach allows us to construct discrete breathers with a 
high level of accuracy for different scalar models. 

We used the same invariant manifolds from Table 1 for constructing DBs for both model 1 
and model 2. For example, for the manifold Q^ 5 , we obtained the following initial conditions 
which determine DB with T = 2 for the model 2 with 7 = 1 (all initial velocities of the 
dynamical variables are assumed to be zero): 

a(0) = 2.23751, 6(0) = -0.64719, c(0) = 0.30201, d(0) = 0.13725, 
e(0) = -0.07879, /(0) = 0.02294. 
On the other hand, using the above discussed method for the model 1 (7 = 1, m = 4), 
we found the following DB with T = 2 on the same manifold Q^5 : 

a(0) = 1.28198, 6(0) = -0.28941, c(0) = 0.00366, d(0) = 0.00180, 

(46) 

e(0) = -0.00018, /(0) ~ 10" 14 . 
(This breather was also obtained with the aid of the Rosenberg nonlinear normal modes 
technique discussed in Sec. IIIII) . 

Comparing solutions (145]) and (|46l) . we see that discrete breather in the homogeneous 
potential model 1 turns out to be much more localized in space than that in the model 2. 

Because of this reason, DBs on the invariant manifolds Q^5 an d Q^x5 f° r the m °del 
2 are essentially different, while DBs on the same manifolds for the model 1 prove to be 
practically identical. Indeed, the peripheral dynamical variables are very small in this case 
and the breathers on and Q^ 5 manifolds are practically indistinguishable from each 
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other up to the used accuracy (note that Q^ 5 transforms into Q^s if g = f = e). 

Surprisingly, using the above described approach, we were able to find some "unusual" 
DBs both for the models 1 and 2. As an example, in Fig. [61 we present DB with T = 2 on 
the manifold Q^ 5 for the model 1 with homogeneous potential of m = 4 degree. Note that 
this discrete breather is not a Rosenberg nonlinear normal mode. Indeed, it can be seen in 
Fig. O that different particles pass their equilibrium position at different instants and this 
contradicts to the definition of the Rosenberg NNMs [see Eq. ([15]) ]. As to our knowledge, all 
papers devoted to discrete breathers in the models with homogeneous potentials deal only 
with the solutions constructing on the basis of separation of space and time variables. Such 
solutions represent Rosenberg nonlinear normal modes. However, one can see in Fig. [6] that 
in the models with homogeneous potentials there can also exist essentially another types of 
discrete breathers. 

The unusual DB presented in Fig. [6] is determined by the following initial values of 
dynamical variables: 



a(0) = 1.6208, 6(0) = -1.3718, c(0) = 0.6325, d(0) = -0.0370, /(0) = -0.0354, 




FIG. 6: Time evolution of DB for the model 1 on the invariant manifold Q5J5. This breather is 
not a Rosenberg nonlinear normal mode. 

while the breather depicted in Fig. [7] is generated by Eq. (1461) . 

In both cases, initial velocities of all the particles are equal to zero. 
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FIG. 7: Time evolution of DB for the model 1 on the invariant manifold Q^b- This breather is a 
Rosenberg nonlinear normal mode. 

Moreover, we can also construct DBs of the above discussed unusual type for the model 
2, as well as for other scalar models on square lattice. 

V. GENERAL GROUP-THEORETICAL METHOD FOR SIMPLIFYING STABIL- 
ITY ANALYSIS OF PERIODIC AND QUASIPERIODIC REGIMES 

A. General discussion 

The system ( |32|) of variational equations for the model 1 possesses a very specific structure 
which permits one to decompose it into independent scalar equations by diagonalyzing the 
matrix B. Unfortunately, the variational equations for studying stability of breathers in the 
model 2 possess more complicated structure which does not allow such decoupling. Indeed, 
the Jacoby matrix J(t) from the variational equations 6 = J(t)S for the model 2 can be 
written as follows 

J(t) = D(t) + C. 

Here D(t) is a certain diagonal matrix whose elements are different time-periodic func- 
tions determined by the breather solution, while C is a constant matrix depending on the 
parameter 7 which describes coupling between Duffing oscillators. 

It is essential, that these two matrices, D(t) and C, don't commute with each other: 

D{t) -C ^ C-D(t). 
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Therefore, it is impossible to diagonalize them simultaneously, i.e. with the aid of one and 
the same orthogonal transformation. 

In turn, this means that one cannot decompose the system of variational equations 
8 = J(t)8 into independent scalar equations. 

Then, one may ask: "Are there any methods for decomposing the system of variational 
equations into some independent subsystems of less dimension than that of the whole system 
5 = J(t)6"? A general group-theoretical method was developed for solving this problem in 



our paper 



341 ] . Below, we give an outline of this method adapted to the 2D scalar models. 



Let us suppose that the system of nonlinear dynamical equations (j7J) is invariant with 
respect to a group Gq. This means that after the action of induced operators g, associated 
with each g G G , we obtain a new system of differential equations which is equivalent to 
the original system (J7J). 

It can be seen from the structure of Eqs. © and the group G (j2J) that under the action 
of any operator g (g € Go) the individual equations are transposed and change their sings 
exactly in such a way, as the corresponding dynamical variables qij. 

On the other hand, as was discussed in Sec. ??, every dynamical regime in the considered 
model can be associated with a certain subgroup G of the parent symmetry group Go 
(G C Go). For studying its linear stability, we linearize nonlinear equations (J7|) near this 
regime and obtain a variational system 8 = J(t)8, where 8 = 8(t) is an infinitesimal vector 
describing perturbations of all the variables qij(t), while J(t) is the Jacoby matrix of the 
system ([7]). 

Then the following question arises: "What one can say about symmetry group of the 
linearized (variational) system 8 = J(t)8?" . It was proved in 34J that this system is 
invariant with respect to the group G describing the symmetry of the given dynamical 
regime. 

Since each operator g generates a transposition (and, may be, changing in sign) of indi- 
vidual variables q^, one can associate a certain matrix with every g £ G. In this way, we can 
construct a matrix representation T of the group G. In general, this representation turns 
out to be reducible and, therefore, it can be decomposed into a number of irreducible rep- 
resentations (irreps) Tj of the group G with the aid of standard group-theoretical methods 



(see, for example, 



33|): 
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(47) 



In 



34J, we proved that the Jacoby matrix J(t) commutes with all matrices of the rep 



resentation T. This fact allows one to apply then the well-known Wigner theorem 33 ]. 
As a consequence of this theorem, the matrix J(t) can be represented in the following 
block- diagonal form 

J(t) = J2 (B D J , (48) 
j 

where each block Dj corresponds to the irrep Tj of the group G. Dimension of Dj is equal 
to rrij ■ rij, where rij is the dimension of the irrep Tj, while rrij determines how many times 
Tj enters into the decomposition ( 1471) of the reducible representation T. 

Moreover, each block Dj possesses a very specific structure: it consists of subblocks 
representing matrices proportional to the identity matrix I n . of dimension rij which are 
repeated rrij times along the rows and columns of the block Dj. 

We can illustrate the structure of a certain block Dj = D characterized by the numbers 



rij = n, rrij 



m as follows 



D 



^2\In l^22^n ■■■ ^2m^n 



(49) 



where I n is the n x n identity matrix. 

As was already noted, individual nonlinear equations ([7]), as well as those of linearized 
system 8 = J(t)8 are transformed under the action of operators g (g G G) as the cor- 
responding variables g^. Therefore, we can construct the same reducible representation T 
of the group G using as the basis the set of equations entering the system S = J(t)S. In 
turn, this leads to representing this system in the form of independent subsystems of rrij ■ rij 
equations corresponding to the individual irreps Tj. Moreover, the specific structure ( )49l) 



of the blocks D~ in (1451 allows one to conclude (see details in 



34j) that the subsystem 



corresponding to the block Dj can be splitted into rij new independent identical subsystems 
whose dimension is equal to rrij. 

Thus, as a result of the decomposition of the variational system 8 = J(t)8, we can obtain 
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for each irrep Tj, entering the representation T, rij identical subsystems whose dimension 
is equal to rrij. 

Obviously, such a decomposition may be very effective, if the integer numbers rrij are 
much smaller than the full dimension of the system 6 = J(t)8. 

It is very important that the above described decomposition does not depend on the 
character of the considered dynamical regime — it can be periodic, quasiperiodic or even 
chaotic! 

As was discussed in [34(, for decomposing (splitting) the original system S = J(t)S into 



independent subsystems in explicit form, we must find the basis vectors of all Tj entering T. 
These vectors can be used as columns of the matrix S which allows one to split the system 
6 = J(t)S. 

Using this matrix, we can introduce a new infinitesimal vector y instead of the old vector 
6 with the aid of relation 

y = SS. (50) 

Note that S is an orthogonal matrix: S T ■ S = I (S T is the transpose of the matrix S). 
Multiplying both sides of the equation 8 = J(t)S by the matrix S from the left, we obtain 
S5 = (SJ(t)S T ) ■ (SS), or 

y = J(t)-y. (51) 

Here the matrix J(t) = S ■ J(t)S T possesses the above discussed block-diagonal structure 
and, therefore, Eqs. (I5TI) represent a set of independent subsystems whose dimensions are 
determined by rrij — the number of times that the irrep Tj enters the reducible representa- 
tion r. 



B. Splitting schemes 

The explicit decomposition of the linearized system d = J(t)S represents a cumbersome 
procedure and, therefore, it is interesting to know beforehand to what extent such decompo- 



sition will be useful. For this purpose, we construct the so-called splitting schemes |34| with 
the aid of the theory of characters of group representations. Let us consider this approach 
in more detail. 

The splitting scheme determines how many independent subsystems one can obtain as a 
result of the decomposition of the linearized system S = J(t)6 and what dimension possesses 
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each of these subsystems. 

To obtain the splitting scheme of variational equations corresponding to a given dynamical 
regime, we must find how many times {rrij) each irrep Tj enters the reducible representation 
T associated with the symmetry group G of the considered regime. 



According to the theory of group representations [33] 



(53) 



Here ||G|| is the number of elements of the group G, Xr(g) and Xj(g) are traces of 
the matrices associated with the group element g G G in the reducible representation Y 
and in the irreducible representation Tj, respectively [the bar over Xj(g) denotes complex 
conjugation]. 

As an example, let us consider constructing of the splitting scheme for a dynamical regime 
on the invariant manifold Q^x3 corresponding to the point symmetry group C± v = {C4, a y }. 
This group consists of the following eight elements (\\G\\ = 8): 

g%(x, y) — (x, y) — identity element; 
g 2 (x, y) = (—y, x) — rotation by 90°; 
gs(x, y) = (— x, —y) — rotation by 180°; 
g4,(x, y) = (y, —x) — rotation by 270°; 

g$(x, y) = (—x, y) — reflection in the mirror plane (1, 0) [a x ]; 

ge(x, y) = (—y, —x) — reflection in the mirror plane (-1, 1) [<r i2/ ]; 

g 7 (x, y) = (x, —y) — reflection in the mirror plane (0, 1) [<J y \, 

g%(x, y) = (y, x) — reflection in the mirror plane (1, 1) [cr xy \; 
Here we define every symmetry elements by its action on an arbitrary point (x, y) of the 
two-dimensional plane. All rotations are performed about Z axis perpendicular to this plane 
and passing through the origin of the coordinate system. Mirror planes are determined by 
their normals which are given as two-dimensional vectors at the end of the corresponding 
lines. (In square brackets, we give notations of these planes used in the previous sections of 
the present paper). 

Let us illustrate definitions of the symmetry elements from Eqs. ( 1531) with the example of 
g 2 . This element, being the rotation by the angle 90°, maps the point (x, y) onto the point 
(— y, x), i.e. x — > —y, y — > x under the action of g 2 . 
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Every dynamical regime on the invariant manifold Q^ 3 is described by three variables 
a(t), b(t), c(t). However, an arbitrary perturbation of this regime is characterized, obviously, 



by nine independent variables 5i 



1..3, j = 1..3): 



Su §12 C>13 ' 



°21 ^22 ^23 
031 °32 033 



(54) 



Thus, we must consider the nine- dimensional space of all possible perturbations 
In this space, we choose a natural basis {ex, e 2 , ...,e 9 }: 

^ l\ 




^10 0^ 



ei 








^0 10^ 



e 2 








e 3 







^0 ^ 



e 9 



V 




1 



Only one component of each vector Bj possesses nonzero value (it is equal to unity). 

Basis vectors ej (j = 1..9) transform into each other under the action of any symmetry 
element g E G. For example, g 7 e x = e 7 , g 7 e 2 = e 8 , g 7 e 3 = e 9 , g 7 e 9 = e 3 . 

This transposition of the basis vectors determines a 9 x 9 matrix M(g 7 ) of the reducible 
permutational representation T corresponding to the element g 7 . In such manner, we can 
construct the complete representation T as the ordered set of matrices {M(g) \ Wg G G}. 

However, for using Eq. ( 152)) . we need only the character of the representation T, i.e. the 
ordered set of traces xid) °f & U hs matrices M(g). 

It is easy to see that nonzero contribution to x{q) originates only from the basis vectors 
Cj which are not changed under the action of the element g. For g 7 such vectors contain 
1 in those sites of the lattice which lie on this mirror plane: <? 7 e 4 = e 4 , ^765 = e 5 and 
g 7 e 6 = e 6 . Therefore, xi&r) = 3 - 

Continuing in such a way, we obtain traces x(g) f° r all g E G = C 4 „ for the reducible 
representation T (see the last line in Table 2). 

Table 2. 

Irreducible representations Tj (j = 1..5) of the group G = G± v and the character x(H 
of its permutational representation V. 
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To 



91 




92 



1 
-1 



03 



-1 
-1 



94 



-1 

-1 

1 



05 



-1 

-1 
1 



9e> 



-1 




97 

1 

-1 
1 

-1 
1 
-1 



9s 



1 

-1 
-1 



x(r 5 



o 



o 



















x(r) 



9 



1 



On the other hand, the group G — C^v 

possesses five irreducible representations: four 
of them are one-dimensional (ri, I^, 1^, r^, while Ts turns out to be two-dimensional. 
This information is presented in Table 2 and can be found practically in all textbooks on 
the group theory. 

Using Table 2 and Eq. ( 1321 . we can easily obtain: 



m(ri) = 3, m(r 2 ) = o, m(r 3 ) = i, m(r 4 ) = i, m(r 5 ) = 2. 

Therefore, the reducible representation T can be decomposed into irreducible representations 
as follows: 

r = 3Vt © r 3 © r 4 © 2r 5 . (55) 

According to the above-discussed method for decoupling the linearized system 8 = J(t)S 
into some independent subsystems, we, therefore, find the following splitting scheme: 



(3 x 3)©(1 x 1) © (1 x 1) ©2(2 x 2). 



(56) 



Thus, 9x9 system of variational equations d = J(t)S can be decomposed into five 
independent subsystems: 

- one of them is three-dimensional (it corresponds to the irrep Fx), 

- two represent different scalar equations, i.e. lxl system, corresponding to the irreps 
T 3 and T 4 , 

- two subsystems are two-dimensional with the same 2x2 matrix (they correspond to 
the irrep r 5 ). 

Obviously, the study of stability of the considered dynamical regime becomes considerably 
easer: instead of analyzing 9x9 system of variational equations one can study stability 
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properties of the four independent subsystems whose dimensions are equal to 1, 2 and 
3. Such approach demonstrates its particular effectiveness in the cases where dynamical 
regimes are quasiperiodic and, therefore, we cannot use the Floquet method for studying 
their stability. 



C. Explicit decomposition of the variational equations 

Now we will consider the method for obtaining the explicit form of the decoupled sub- 
systems of the linearized system 6 = J(t)6. 

As was already discussed, to this end, we must find the orthogonal matrix whose columns 
are basis vectors of all the irreps Vj entering the reducible representation T. Usually, basis 



vectors of irreducible representations of a given group are constructed 



the projection 



operator method [33|. However, we prefer to use the "direct" method [31] based on the 
straightforward application of the group representation definition. This method turns out 
to be more obvious and more simple for our purposes. 

Suppose $ = {Vi, V2, V n } is a " supervector" representing the set of all basis vectors 
Vj of a given (reducible or irreducible) n- dimensional representation V of the group G: 

T = {M(g)\\JgeG}. 



This means that (see, for example, 33(): 



g<f> = M T (g)<f>, VgeG, (57) 

where M T (g) is the transpose of the matrix M(g) associated with the element g of the group 
G. 

We want to construct basis vectors of all irreps Tj (j = 1..5) of the group G = C^ v 
presented in Table 2. Since 1^, T 2 , T 3 and T 4 are one- dimensional, only one basis vector 
corresponds to each of them. We denote these vectors by <pi, cp 2 , ¥>3 and cp 4 , respectively. 

Note that it is sufficient to take into account only those equations from fl57|) which cor- 
respond to generators of the group G: all other elements can be constructed as a different 
products of these generators and equations (I5T1) . written for such elements g G G, don't give 
us any additional information about the basis vectors of a given representation. We choose 
elements g 2 and g^ as generators of the considered group G — C^ v . 
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Then, from the definition (|57|) . we obtain the following equations for the above basis 
vectors: 

Ti : g 2 tpx = Vii faVi = Vi\ (58) 



r 



2 : 92^2 = <P2, 95<f2 



-<P2\ 



(59) 



3 : g2fi 



"V3, 95Pl = <f3; 



(60) 



92<f4 



-<P4, 95<P4 



(61) 



On the other hand, the irrep T 5 is two-dimensional and, therefore, two basis vectors, Vi 
and ?/>2, correspond to it. They satisfy the following equations: 



(62) 



g 2 1pl = -02, 92^2 = -V>i! 
05^1 = -il>i, 9r,4>2 = ip2- 

Solving Eqs. ([58]) - ([52]) . we can construct the complete set of basis vectors in nine- 
dimensional space of all possible deviations %(£) [see Eq. fl54"j) ] from the considered dynam- 
ical regime on the invariant manifold Qsx3- 

Noting that Eqs. fl58|) coincide with Eqs. (fl~2"]) - ffT3]) up to renaming variables (g^ — > 5ij), 
we conclude that the basis vector of the identical irrep Ti possesses the same form as the 
invariant manifold Q^3 i n Eq. ffT4")) : 



f C B 



BAB 
C B C 



(63) 



Thus, (fi depends on three arbitrary parameters A, B, C and represents a three- 
dimensional subspace of the above-mentioned nine- dimensional space. This means that 
the irrep Ti enters the reducible representation V three times. 

Taking into account Eq. fl63|) . we can rewrite ipi as follows: 



/ 



<Pi = A 




1 




\ 



B 



1 

1 1 
1 



\ 



/ 



c 



1 1 \ 

10 1/ 



(64) 
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Therefore, one can choose the following three basis vectors Vj in the three-dimensional 
subspace determined by ifi. 



Vi = (0,0,0|0,1,0|0,0,0) T ; 
F 2 = ±(0,1,0|1,0,1|0,1,0) T ; 

y 3 = |(i,o,i|o,o,o|i,o,i) T 



(65) 



These nine-dimensional vectors form the first three columns of the matrix S which leads to 
the decomposition of the linearized system 6 = J(t) ■ S into some independent subsystems. 
(Here we use the transposition symbol T to present column vectors as row vectors). 

From Eqs. (I59p . we find if2 = 0, i.e. the irrep T 2 does not enter the reducible represen- 
tation r. 

The solution 

/ -B 
BOB 
0-5 

of Eqs. (1601) means that the irrep T 3 enters the representation T only one time and, therefore, 
the fourth column of the matrix 5" represents the vector 

V 4 = i(0,-l,0|l,0,l|0,-l,0) T 



<P3 



V 



\ 



(66) 



From the solution of Eq. f l6"T|) 



(f 4 



( -A A ^ 


A -A 



(67) 



we obtain the fifth column of the matrix S: 



1 



V 5 = -(-l,0,l|0,0,0|l,0,-l) J 



Finally, for the irrep T 5 , we find from Eq. 

1 A B A ^ 



V 



ABA 

-A -B -A 



V>2 



( A -A ^ 
B —B 
A -A 



(68) 
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These basis vectors depend on two parameters (A, B) and, this confirms that the irrep Y§ 
enters the representation Y two times [see Eq. fl55|) ]. Supposing successively A = 1, B = 
and A = 0, B = 1, we obtain after normalizing the last four columns of the matrix S: 

F 6 = ±(1,0,1|0,0,0|-1,0,-1) T ; 
V 7 = ^(0,1,0|0,0,0|0,-1,0) T ; 
V 8 = ±(1,0,-1|0,0,0|1,0,-1) T ; 
F 9 = ^(0,0,0|1,0,-1|0,0,0) T 

All the vectors Vj (j = 1..9) are orthogonal and normalized. Note that this is the 
necessary condition for the matrix S to be orthogonal. 

As a result of all above-discussed steps, we confirm the correctness of the decomposition 
( )55l) of the reducible representation T into the irreps of the group G = C± v and construct the 
following matrix S which allow us to split the linearized system 8 = J(t) ■ d into independent 
subsystems: 



S 







i -i 
\ 

± o h 



2 

1 

o I 



I 

o I o -I 



1 

2 






1 

2 



1 

2 












1 

2 



^ 











- 



v/2 
2 



i -i -i -i 







yj 

2 



2 



\ 









/ 



(69) 



As was described in Sec. IV Al [see, Eqs. ( l50l) -(l5TT)]. introducing the new nine-dimensional 
infinitesimal vector y = S8 with matrix S determined by Eq. ( |69|) . we obtain from the old 
linearized system 6 = J{i) ■ 5 the new one 



y = J(t)y, J(t) = SJ(t)S 7 



(70) 



which turns out to be decomposed into some independent subsystems. 
Example 1 

For the model ([2]) of linearly coupled Duffing oscillators, we obtain the following subsys- 
tems corresponding to the individual irreps of the group G — C^ v . 
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The irrep T± generates the following 3x3 subsystem: 



Vi + [3c 2 (t) + (1 + 2a)\y x - 2ay 2 = 0, 

y 2 - 2ay x + [3b 2 (t) + (1 + 3a)]y 2 - 2ay 3 = 0, 

y 3 - 2ay 2 + [3a 2 (t) + (1 + Aa)]y 3 = 0. 



(71) 



All three time-periodic functions, a(t), b(t) and c(t), describing a dynamical regime on 



' c b c ^ 



the invariant manifold Q 



(i) 

3x3 



V 



6 a b 
c b c 

that this regime can be periodic, as well as quasiperiodic. 
Two identical 2x2 subsystems 



enter Eqs. (|71|) . Let us emphasize once more 
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y x + [3c 2 (t) + (1 + Aa)} Vl - V2ay 2 = 0, 
y 2 - v^ayi + [3& 2 (t) + (1 + 5a)]y 2 = 



(72) 



correspond to the irrep T^, while different scalar equations correspond to the irreps 1^ and 
T 4 , namely, 

r 3 : y+[3b 2 (t) + (l + 3a)}y = 0, (73) 



T 4 : y+[3c 2 (t) + (l + 6a)]y = 0. (74) 

For studying stability of discrete breathers in the model 2, we can apply the standard 
Floquet method. In the case of the invariant manifold Q^ 3 , one has to deal with the 
18-dimensional phase space corresponding to a given breather. On the other hand, after 
splitting the linearized system <5 = J(t) ■ 5 into independent subsystems (|7T]) - fl71|) . we 
may analyze discrete breather stability applying the Floquet method successively for these 
subsystems of much smaller dimensions. 

For DB with T = 2 determined on the invariant manifold Q 3 ^ 3 in the case of the model 
2 with 7 = 1 by the initial conditions, a(0) = 2.35593, 6(0) = -0.45156, c(0) = 0.13152 
we have obtained the following Floquet exponents for 3 x 3 subsystem (ITT!) : 

-0.74678 ±0.66507 ■ i, 
< 0.26213 ± 0.99966 ■ i, (75) 
l±0-i. 
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For 2x2 subsystems (172 



-0.60978 ± 0.79257 ■ i, 

(76) 

0.61449 ± 0.78892 ■ i. 



- 0.53267 ± 0.84632 ■ i. (77) 



For Eq. fl73]): 
For Eq. (1741): 

- 0.55546 ± 0.83154 ■%. (78) 

All these Floquet exponents lie on the unit circle in the complex plane and, therefore, 
the given DB proves to be stable. 

To check the above-discussed splitting procedure, we calculated the Floquet exponents 
for the whole 18-dimensional phase space associated with the original linearized system 
5 = J(t)S and are convinced that they exactly coincide with those from Eqs. (1751) - (1781 . 

Actually Eqs. (177)1) - (1781 allow us to study stability of a given DB with respect to different 
symmetry-determined collective degrees of freedom. In some a sense, such information can 
throw light upon the cause of the breather stability loss. 

Example 2 

Finally, let us decompose the linearized system d = J(t)d for the discrete breather on 
Q3X3 in the case of the model 1 whose stability we were able to study without applying 
the Floquet method (see Sec. IHIj) . Using the same matrix S from Eq. (1691) we obtain the 
following decomposition of the linearized system: 

m + 3[ 7 c 2 (t) + 2u(t)] yi - 6u{t)y 2 = 0, 
r, : {y 2 - 6u(t) Vl + 3[jb 2 (t) + M (t) + 2u(t))y 2 - 6^(t)y 3 = 0, (79) 
y 3 - 6fi{t)y 2 + 3[ 7 a 2 (t) + 4/i(£)]y 3 = 0, 



5 • 



y x + 3[ 7 c 2 (t) + 2u{t)} Vl - 3V2u{t)y 2 = 0, 
y 2 - 3V2v(t) yi + 3[ T 6 2 (t) + fi(t) + 2v{t)\y 2 = 0, 
(two identical systems of this form), 



(80) 



y + 3{ 1 b*(t) + f ,(t) + 2v(t)]y = 0, (81) 
: y + 3[ 7 c 2 (t)+2iy(t)]y = 0. (82) 
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Here fi(t) = [a(t) - 6(f)] 2 , v(t) = [&(*) - c(t)] 2 , where a(t), b(t), and c(t) are three breather 
components on the invariant manifold Q^3- 

As a consequence of this decomposition, we can associate different stability indicators Aj 
(j = 0..8) (see Sec. IIIip for DB (7 = 1) with the irreducible representations r 1; T 3 , T 4 , T 5 
of the group G = C± v : 





Ao 


= 3, Ai = 0.12407, 


A 2 = 0.02340 


r 3 : 


A3 


= 0.59433 




r 4 : 


A 4 


= 0.03707 




r 5 : 


A 5 


= A 6 = 0.59556, A 7 


= A 8 = 0.03584 



(83) 



Note that some of these Aj turn out to be degenerate, exactly or approximately. We will 
consider symmetry-related causes of this phenomenon elsewhere. 

In connection with the above group-theoretical method for splitting the linearized sys- 
tem 8 = J(t)S into independent subsystems, one may ask: "How this method works in 
the case of large lattice fragment?" Indeed, the worse the localization of the breather, the 
larger fragments of the corresponding invariant manifolds one must consider. Unfortunately, 
our decomposition method becomes less effective for this case because dimensions of the in- 
dependent subsystems can be rather large. For example, we obtain the following splitting 
scheme for studying stability of DB with symmetry group G = C^ v on the invariant manifold 
{1) ■ 

1(6 x 6) © 1(1 x 1) © 1(3 x 3) © 1(3 x 3) © 2(6 x 6). (84) 

These subsystems are associated with the irreps Ti, T 2 , r 3 , T 4 and T 5 , respectively. 

Comparing ( 18~4"]) and (15"6"1) . we see that the decomposition of the linearized system on the 
manifold is more efficient than that on Q^ 5 . 

However, let us emphasize again that all discussed decompositions are true not only for 
periodic, but for any quasiperiodic regimes, as well. In practice, the latter can be very 
important. 

VI. CONCLUSION 

In the present paper, we study stationary discrete breathers in two nonlinear scalar 
models on the plane square lattice. However, the developed group-theoretical methods for 
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constructing breather solutions and studying their stability can be naturally extended to 
the case of time-periodic and quasiperiodic nonlinear dynamical objects of various physical 
nature on different 2D and 3D lattices. These methods allow one to simplify, sometimes 
considerably, studying discrete breathers and quasibreathers in space-periodic structures. 
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